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pLn ' Color image segmentation is an important topic in the image processing field. MRF-MAP is often 

I adopted in the unsupervised segmentation methods, but their performance are far behind recent interactive 

segmentation tools supervised by user inputs. Furthermore, the existing related unsupervised methods 
■ also suffer from the low efficiency, and high risk of being trapped in the local optima, because MRF-MAP 

is currently solved by iterative frameworks with inaccurate initial color distribution models. To address 
Q , these problems, the letter designs an efficient method to calculate the energy functions approximately 

in the non-iteration style, and proposes a new binary segmentation algorithm based on the slightly 
tuned Lanczos eigensolver The experiments demonstrate that the r 
performance compared with two state-of-art segmentation methods. 



^ ■ tuned Lanczos eigensolver The experiments demonstrate that the new algorithm achieves competitive 
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■ Image segmentation, Markov random fields, maximum a posteriori, unsupervised segmentation. 

^ . I. Introduction 

H. 

. Unsupervised color image segmentation is important in various image processing and computer vison 

applications, such as medical imaging HI, image retrieval 121, image editing Q, and object recognition 
im. Estimating the maximum a posteriori (MAP) on the Markov random fields (MRF), is so far an 
fundamental tool which is widely adopted both in the unsupervised color image segmentations and 
supervised ones |[5l- |[T0l . In all existing MRF-MAP-based image segmentation methods, their goals are 
to find the optimal label configurations on pixels to maximize the posterior probabihty which is in 
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proportion to the product of the MRF priors and the likelyhooods terms, or equivalently, to minimize the 
energy function of the smoothness terms plus data terms. 

The Ukelyhood/data terms consist of the parameters of some color distribution models, here these 
models specify the probabilities of any color occurring in each segmentation. These parameters usually 
derive from user interactions or random sampling, hereby there is always non-negligible inaccuracy in 
the lUcelyhood terms. In order to address this, the Expectation-Maximization (EM), simulated annealing 
and other iterative methods are usually adopted to progressively approach the appropriate parameters, 
especially in the case of unsupervised segmentations up to today m-Q. There are many choices of 
optimization algorithms to be adopted in the M-step ifTTl . 

There are three major disadvantages in the current unsupervised segmentation methods. The first is the 
low efficiency of the iterative frameworks, particularly when faced with large size images. The second is 
always on the high risk to be trapped to the local optima of the energy function. Although stepping out and 
restarting the iteration is a reasonable improvement, there would inevitably be additional computational 
load and it is possible to be trapped again. The third is the coarseness in the segmentation results, and 
it is partially caused by the roughness of the likelyhood parameters. To address these issues, the letter 
proposes a new unsupervised binary segmentation method based on the approximation of the likelyhood 
terms, where the iterative computations are replaced by only one single step of solving the eigenvector of 
the largest eigenvalue, therefore the computational efficiency is remarkably improved. This new approach 
increases the chance to high quality segmentation results by obtaining the nearly optimal solutions to 
maximize the posterior probability as possible. It also provides us an effective way to test and verify the 
MRF prior parameters or their involved generating schemes, which are also critical to the segmentation 
tasks. 

II. Color Image Segmentation based on MRF-MAP 

In the following sections, we focus on the binary segmentation with the label set {fore, back}. The 
computational goal of these methods is to maximize the probability P{L\I) of the segmentation label 
configuration L given an image /. According to the Bayesian rule, it is equivalent to maximize the joint 
probability P{I, L) = P{I\L) ■ P{L), where the prior P{L) is established on the Markov random field of 
/, the conditional probability P{I\L) is the likelyhood that the pixel colors occur in their corresponding 
segments marked by different labels. In a more prevailing view, what we need is to minimize a energy 
function E which is the negative log-hood of P{I,L). Here E is usually written in the form of a data 
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term Eo plus a smoothness term Es multiplied by a factor A: 

E = Ed + XEs (1) 

where Ed reflects the likelyhood of the color occurrences in the image segments, Es is the sum of all 
adjacency interaction potentials of each two neighboring pixels of different labels: 

Ed = ^ -In P^p){p),Es= J] S(c(p), c(g)) (2) 

P ip,q)eN 
L{p) + L{q) 

where L(-) is the pixel label, c(-) is the pixel color, and S{-, •) is the perceptually similarity weight of 
two colors. It is meant two pixels p and q are adjacent to each other by noting (p, q) € N. Thereafter, the 
segmentation task is to pursue an appropriate a label configuration to reach the lowest energy. Although 
as indispensable as the data terms when computing the energy functions, the smoothness terms are not 
to be addressed in the letter. 

There are two steps when determining the coefficients in the data term Ed- First, a suitable color 
statistical model should be chosen. Histograms are usually adopted for images of small color spaces, 
such as gray scale or 256 colors, but it is not suitable for large color spaces as the samples were 
statistically too few when facing so many histogram bins. Some other models fit large color spaces well, 
and make an appropriate comprise between the efficiency and accuracy, such as the Gaussian mixture 
model (GMM). Second, the model parameters should be determined. However here arises a chicken or 
the egg dilemma unavoidably: we have to know the parameters first to minimize the energy function to 
obtain the optimal segmentation, but the optimal segmentation is just the key to produce the accurate 
parameters mentioned above. The usually adopted solutions to this are the iterated procedures, such as 
EM, in which the estimation and optimization are performed sequently but isolatedly in each single loop. 
Here the initial parameters are determined from sample pixels chosen by user interactions or random 
samplings. There are many choices to perform the optimization ifTTl : graph cut. Loopy Belief Propagation 
(LBP) and Iterated Conditional Model (ICM). 

The low computational efficiency is an adherent shortcoming of the iterated solutions, as it is extremely 
hard to predict when and where the iterations would halt. Furthermore, although the minimum of 
the energy function is an unambiguous target itself, the actual aim of the iterated solutions is not 
mathematically explicit for us to approach. As a result, the iterations are likely finished at the local 
minimal in most cases. To address these issues, the letter proposes an approximating expression which 
is rather close to Ed in (2), and associate the approximated target energy function with the cut on a 
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complete graph G which has both positive and negative edge weights. In the following manipulations, it 
is rather straightforward to solve an eigen-system to pursue the minimum cut C on G, so to minimize the 
energy. There the expected segmentations are worked out directly without considering the troublesome 
parameters of the data terms at all. 

III. Segmentation based on Approximate MRF-MAP 

Consider there are n pixels of m{m«n) colors in the image /, and any two colors are perceptually 
distinguishable from each other. It is almost always achievable with the help of the existing color 
clustering algorithms, even facing much splendid images. Therefore we choose the histograms as our 
color distribution models in (2). Let c(p) be the color of the pixel p. For all pixels having color i, let 
Ui be the total amount, n /^^ be the amount of those having label fore, so be nt^j and label back. Clearly 
rii = nf^i + j. Let F and B be the two sets of fore and back pixels respectively, and the corresponding 
pixel amounts are n/ = nf^i + nf^2 + ■ ■ - + 1^ f,m and ni, = n^^i + nb^2 + • • • + nfe,m- Leaving the smoothness 
term Eg unchanged, we have the the data term Ed as 

Ed = y -In^^ + y -In^^ 

m 



(n/lnn/ + n6lnn6)-^(n/,jlnn/,j + nb,i\n.nb,i) 



i=l 



(3) 



n n n n 



[rii ■ (^^In^^ H ^In— ) + njlnnj 

Tli Tli Tli Tli 



i=l 



Now consider a function g{x) defined on the interval [0, 1]: 

{a;lna; + (1 — x) ln(l — x),if < a; < 1. 
(4) 
0, otherwise. 

Clearly g is continuous on the whole interval [0, 1]. Expanding the log terms into the Taylor series and 
simplifying the expression, we have g{x) = — fxy + A(a;), where y = (1 — x) and 

A(x) = -xy[^{x'^ + y'^) + ^{x^ + y^) + ■■■], (5) 

which has the mean value A{x)dx = — on [0, 1]. If we replace A(x) with the constant — on 
the entire interval, the mean squared error would be (A(a;) — < 3 x 10^^. Since the MSB is 

considerable small, it is totally acceptable to approximate g with g* = —^xy — ^ shown in Fig. 1. 
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Fig. 1. Comparison of function g and its approximation g* on [0, 1]. 



With the help of g*, now we can approximate the energy function E in (3) with 

5 1 " 5 1 

{-—nfnb~ — + nlnn)-'y] {-—nj^i-rib^i ~ — +njlnni) 



1=1 

(6) 



(p, q)<^N 
L(p) + L(q) 

Since m, n and all nj's are constant in the input image, it is actually to minimize the following function 

E* when minimizing the above expression as the approximation of E: 

5 ""5 
E*=-iz-nf-nb+Y.{l^nf^i-ni,^i) + \-'Y^S{c{p),c{q)). 

(p, q)eN 
L{p) + L{q) 

Now construct an undirected complete graph G of n nodes each of which corresponds to one pixel, 
and set the edge weight w{p,q) to be the sum of the following three terms 

2I-, if c{p) = c{q) = i, 

0, otherwise. (8) 



W2{p, q) 
W3(.P, q) 



X-Sic{p),ciq)), if {p,q)€N, 
0, otherwise. 

It is easy to prove that, for any label configuration L of the binary segmentations, E* is equal to the 
capacity of the cut C = {F, B} arising from L on G. Therefore it is equivalent to find the minimum cut 
on G when minimizing E* . However the existing minimum cut algorithms are not suitable here because 
of the existence of negative edge weights. In fact, the problem here is computationally equivalent to a 
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well-known NP-complete problem, Max-Cut, on the graphs with non-negative weights. Therefore it is 
hard to obtain the exact minimum of E* in polynomial time. 

Our solution to this is to generalize it into the continuous real space ii". First we put the label 
configuration into an indicator vector D = [di, ^2, ■ ■ ■ , dn]^'- for i = 1, 2, • • • n, let di = +1 if the ith 
pixel label is fore, and —1 for back. Then establish a matrix W = [w{p, q)] and let Sw denote the sum 
of all its entries. It is easy to prove that the cut value is equal to ^{Sw — D^WD). After generalizing 
the (ij's to be in the continuous interval [—1, 1] instead of {+1, —1}, our task becomes into 

maxD^WD,s.t\\D\\2 = n (9) 

because Sw is a constant here. According to the Lagrange multiplier method |fT2l . the solution to (9) is 
the eigenvector D* corresponding to the largest eigenvalue of W. The Lanczos algorithm, well known as 
the fastest method solving extremal eigenvectors for large sparse matrices, is adopted here to calculate 
D* = [djjfig, • ■ ■ jC^n]^- Since W is full, the embedded matrix- vector multiplication must be improved 
using the special structure of W. At last, we get the required binary labels straightforwardly by setting 
the ith label to be fore if d* > 0, or back if d* < 0. Here is the outline of our new segmentation 
algorithm: 

Algorithm 1 Color image segmentation based on MRF-MAP 
1: Clustering all colors into m classes; compute S{c{p), c{q)) for each pair of adjacent pixels p and q; 

2: Calculate the largest eigenvector D* = [dj, dg, • • ■ , d^f'" of the matrix W with Lanczos eigensolver, 

where we obtain the product R — [ri , r2 , • • • , of W and any vector V = [vi,V2, - ■ ■ , Vn]'^ as 

• let if - 0, and 9i - for alH = 1, 2, • • • , m; 

. for = 1 to n, let 99 ^ - 1^ ■ Vk, 9c{k) ^ Gc(k) + 2^ " '"k, and = Y.j:{j,k)eN ^(j^ ^) ' '^i' 
, for k = 1 to n,rk^ if + 6^^^^) + - ' '"k- 

3: Output the label of the kth pixel as back if d^, < 0, ov fore if d*^ > 0. 



The analysis on the computational complexity is rather straightforward: the step 1, 3 can be finished in 
0{n) time; the matrix-vector multiplication can also be finished in 0{n) time, so the Lanczos algorithm 
revoked in step 2 can be finished in 0(n • d) time, where d is the amount of performed iterations solving 
the eigenvector. Therefore the total time complexity of is 0{n ■ d). Since d is irrelevant with n and 
empirically always less than a certain constant, the new segmentation algorithm is practically a nearly 
linear one. 
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Fig. 2. Comparison of the segmentation results by our method and other two methods. (Leftmost: original images. Other 
columns from left to right: segmentation results by our method with A = 1, 5, 10, by the method in (8) and 1141 .) 

IV. Experiments 

The experiments are finished on color images chosen from two segmentation datasets from Berkeley 
and MS research at Cambridge, together with the source codes developed in Matlab. The color clustering 
method adopted in step 1 is from |fT3l . and all of the amounts of color classes are set to be 16 for these 
images. A varies from 1 to 10, and all smoothness terms are simply set to be 1 over all adjacent pixels 
despite their colors, to say, 5(-, •) = 1. 

These images are intentionally chosen to be of splendid colors and delicate local details, so to verify 
the performance of our new method facing different challenges. Two state-of-art segmentation algorithms, 
one supervised |8| and another unsupervised |[T4l . are chosen in the control experiments to examine the 
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segmentation quality of our new method. In general, all segmentation results of our new method are 
basically acceptable in the experiments. Especially, the perceptually outstanding objects, if any in the 
test images, are usually figured out of the underlaying scenes accurately. There are six groups of typical 
experiment results in Fig. 2: each group includes five different segmentation results of an identical color 
image, three of them are all from our method but with different A's, whereas the rest two are from 
the other two methods as the comparison. The embedded objects, striking either for fresh colors or 
large continuous shapes, are precisely outlined by our new method in Fig. 2. Its segmentation quality, 
is considerably close to the user-interactive-styled method in [8|, but much better than the unsupervised 
one in lfT4l . 

Consequentially, there are some noticeable slight differences between the segmentation results on 
different A's of our new method. There are more isolated, but vivid pieces when A = 1, and simultaneously 
the segmentation boundaries are more likely located on the desired edges of the objects. However, it also 
brings too much emphasis on these discontinuous line segments, and results in much more isolated pieces 
in the segmentations. When A varies from 1 to 5, then to 10, it is shown the segmentation boundaries 
become smoother and smoother on the cost of losing the elaborate details, and the two segmentation 
zones are more close to each other in sizes when A = 10. The reason for that is, the smoothness terms 
become larger and larger quantitatively so that the continuity of the segmentations is emphasized much 
more, therefore it is more inclined to cut two pieces of the the same sizes and flat boundaries to reach 
the minimum energies. 

V. Conclusions 

In this letter, a new unsupervised MRF-MAP-based segmentation algorithm is introduced. By introduc- 
ing an reasonable approximation to the data terms, the energy functions could be minimized remarkably, 
however, without any supervision. The new method is able to obtain the high-quality segmentation results, 
as well as the high computational efficiency. The future work includes the investigating the computational 
hardness of MRF-MAP, and extending the new method to video segmentations and multiple-labeled 
segmentations. 
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